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SECTION  I.  INTRODUCTION 


The  boundary- layer  calculation  method  [1],  which  forms  the  subject 
of  this  verification  study,  was  developed  by  Sybucon  under  Air  Force 
funding  via  a subcontractual  relationship  with  Grumman  Aerospace  Corpor- 
ation. The  method  involved  a radically  different  numerical  scheme  and 
geometry-handling  scheme  than  had  been  employed  in  earlier  methods,  but 
essentially  the  same  turbulence  model  as  one  that  had  been  in  continuous 
use  for  several  years  [2],  This  turbulence  model,  based  on  a still 
earlier  two-dimensional  model  of  Bradshaw  et.  al.  [3],  incorporates  the 
empirically  modified  turbulent  kinetic-energy  equation  together  with  an 
assumed  relationship  between  the  turbulent  intensity  and  the  Reynolds 
shear  stress. 

The  form  of  this  relationship,  and  the  forms  of  two  other  empirical 
functions  of  position  through  the  boundary  layer  (the  dissipation  length 
and  the  diffusion  function),  have  to  satisfy  certain  rather  broad  theore- 
tical constraints  but  otherwise  must  be  adjusted,  by  trial  and  error,  to 
give  the  best  overall  agreement  with  experiment.  A process  of  optimiza- 
tion, on  these  lines,  was  carried  out  six  to  seven  years  ago  [4,  5]  and 
the  resulting  model  has  since  enjoyed  a reputation  for  good  reliability 
and  accuracy.  The  emergence  of  more  recent  sets  of  experimental  data, 
however,  has  indicated  the  potential  for  further  improvements  of  the 
model.  Some  of  these  sets  of  data  include  more  detailed  measurements, 
particularly  of  the  turbulence  quantities,  than  were  available  when  the 
earlier  validation  was  performed.  The  need  for  validation  of  method  in 
compressible  flow  has  also  been  recognized  for  some  time. 
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The  study  reported  here  was  far  from  being  extensive;  comparisons 
were  made  with  only  two  sets  of  measurements  in  three-dimensional  flows. 
There  is  a possibility,  on  that  account,  that  the  proposed  modifications 
to  the  turbulence  model  will  be  biassed.  In  an  attempt  to  minimize  such 
a bias  some  supporting  comparisons  have  been  made  with  a few  sets  of  two- 
dimensional  data.  The  results  are  fairly  convincing,  but  caution  must  be 
exercised  in  interpreting  their  generality.  Until  more  complete  verifi- 
cation has  been  undertaken  the  proposed  changes  must  be  regarded  as 
tentative. 

The  strategy  adopted  for  verification  and  fine-tuning  of  the  turbulence 
model  consisted  of  three  steps: 

“ comparison  with  experiment,  using  the  existing  model 

o determination  of  the  sensitivity  of  the  predictions 
to  changes  of  the  empirical  functions 

° comparison  with  experiment  using  various  candidate 
revised  models. 

Section  II  reviews  the  various  relationships  which  define  the  existing 
turbulence  model,  and  Section  III  presents  the  comparisons  between  this 
model  and  the  experimental  test  cases.  The  sensitivity  study  and  the  pro- 
posed modifications  to  the  model  are  described  in  Section  IV. 
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SECTION  II.  DESCRIPTION  OF  TURBULENCE  MODEL 


The  calculation  method  involves  the  simultaneous  solution  of 
five  governing  equations,  together  with  various  peripheral  relation- 
ships : 

o two  mean-flow  momentum  equations 

o thermal  energy  equation 

o continuity  equation 

o empirically  modified  turbulent 

kinetic-energy  equation. 

This  last  equation,  which  provides  closure  of  the  system  of  mean- 
flow  equations,  is  written  in  the  form 


Pil  .L_ 

hj  8x 


(pV 


3 q^ 

P’v)  37  ( 2 ) 


^ hj  3z  ^ 2 ■* 


where 


puv  3U  ^ pvw  ^ + P_  3_  fP'v  + q^V) 
h2  3y  h2  3y  h2  3y  P 2 


+ pe  = 0, 


“I 

q^  = u'^  + 


~2  T 

V + w 


(1) 

(2) 


Additional  empirical  relationships  are  provided  to  model  the  diffusion 
and  dissipation  terms,  and  to  model  the  dependence  of  the  turbulent 
shear  stresses  and  turbulent  heat  flux  on  q . 

Specifically,  following  Bradshaw  et  al . [3]  and  Nash  [2],  the 
diffusion  and  dissipation  terms  are  modeled,  respectively  by 


Eli 

P 


+ 


e = 


L_ 


(3) 

(4) 
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where  is  the  maximum  value  of  in  the  outer  three-fourths  of  the 

boundary- layer,  and  where  the  diffusion  function:  a^,  and  the  dissipation 
length:  Lp,  are  assumed  to  be  universal  functions  of  physical  distance 
through  the  boundary  layer.  The  existing  forms  for  a.^  and  Lp,  are  the 
same  as  those  found  appropriate  in  earlier  studies  [4,  5],  and  can  be 
represented  by 

= 1.125  - 0.375  n**  (5) 


5 = 7.195n 

6 1 + 4n2  + 5n’ 


(6) 


These  two  functions  are  illustrated  in  Figure  1.  In  Equations  (2.5,  6) 
n is  the  dimensionless  physical  distance  through  the  boundary  layer: 


and  6 is  the  boundary- layer  thickness,  defined  as  the  distance  from  the 
surface  where 


{(U  - U)2  + (W 

e e 


0.005  Q . 


(8) 


Close  to  the  surface  the  dissipation  length  is  modified  by  the  Van  Driest 
correction  [6] 

(Lp)  = L (9) 

corr. 

where 

y*  = f (10) 

and  A is  a constant  usually  taken  to  be  26. 
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The  turbulent  shear-stress  components  are  related  to  by  empirical 


functions  of  the  form 

(11) 

(12) 

In  the  method  of  Reference  1 the  function  f was  tentatively  assumed  to  be 
linear 

f(t^)  = 0.0225  iji,  (13) 


but  the  intention  was  to  study  the  implications  of  using  other  forms  at  a 
later  date.  Some  of  these  implications  are  explored  below. 

The  corresponding  relationship  between  turbulent  heat  flux  and 
is  written  as 


in  which  Pr^  is  an  assumed  constant  value  (=  0.9)  of  the  turbulent  Prandtl 
number . 
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SECTION  III.  COMPARISONS  WITH  EXPERIMENT 


Experimental  Test  Cases 

The  sets  of  experimental  data,  with  which  the  predictions  were 
compared,  are  described  below. 

1 . Two-dimensional  incompressible  flow  over  a flat  plate. 

The  classical  experiment  of  Klebanoff  [7]  provided  the 
main  source  of  data  for  this  test  case.  Mean- flow  data 
and  turbulence  measurements  are  reported  in  Klebanoff s 
paper.  The  mean-flow  measurements  of  Smith  and  Walker 
[8]  and  empirical  correlations  proposed  by  Rotta  [9] 
and  Coles  [10],  provided  additional  data. 

2.  Two-dimensional  incompressible  flow  in  an  adverse  pressure 
gradient . 

Two  sets  of  data  were  selected  to  represent  this  type  of 
flow:  the  measurements  of  Ludwieg  and  Tillman  [11]  in  a 
strong  adverse  pressure  gradient,  and  the  measurements  of 
Schubauer  and  Spangenberg,  "Flow  A"  [12] . These  two  sets 
of  data  were  designated  Numbers  1200  and  4400,  respectively, 
in  Stanford  Conference  proceedings  [13].  Mean-flow  measure- 
ments, only  are  given  in  these  references. 
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3.  Three-dimensional  incompressible  flow  over  a swept  wing. 


The  mean-flow  data  of  van  den  Berg  and  Elsenaar  [14],  and  the 
turbulence  measurements  of  Elsenaar  and  Boelsma  [15]  (which 
were  made  on  the  same  wing  model)  formed  the  input  to  this 
principal  test  case.  The  measurements,  which  were  made  in 
the  Netherlands,  were  designed  to  closely  simulate  infinite- 
swept  wing  conditions.  The  small  residual  variations  in  the 
spanwise  direction  were  accounted  for  in  the  calculations. 
Interesting  background  information,  pertaining  to  the  present 
study  is  provided  by  the  attempts  of  other  investigators  to 
match  the  van  den  Berg-Elsenaar  data;  these  attempts  are 
reported  in  Reference  16. 

4.  Three-dimensional  supersonic  flow  through  a swept  shock. 

The  data  of  Peake  [17]  were  used  for  this  second  principal 
test  case.  The  experimental  arrangement  consisted  of  a 
shock  generator  standing  normal  to  a flat  plate.  The  shock 
formed  by  the  generator  interacted  with  the  boundary  layer 
developing  on  the  plate.  The  measurements  were  made  in 
this  region  of  interaction  where  the  flow  suffered  a strong 
lateral  perturbation.  Two  sets  of  data  are  reported  in 
Peake's  paper;  the  one  at  a Mach  number  of  2 was  selected 
for  comparison  purposes  because  boundary- layer  separation 
was  not  provoked  by  the  interaction.  Separation  did  occur 
in  the  other  set  of  data,  at  M = 4. 
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Comparison  Between  the  Original  Model  and  Experiment 


The  comparisons  between  predictions,  using  the  original  turbulence  model, 
and  experiment  are  presented  in  Figures  2 through  10.  For  the  two-dimen- 
sional flat-plate  flow  (Figure  2)  the  calculation  is  in  good  agreement 
with  the  measured  velocity  and  turbulent-kinetic-energy  profiles.  However, 
the  wall  shear  stress  is  underpredicted:  is  predicted  to  be  0.00256, 

compared  to  the  measured  value  of  0.0028.  The  discrepancy  is  associated 
with  the  failure  of  the  method  to  recover  the  precisely  appropriate  addi- 
tive constant  in  the  logarithmic  wall  law;  it  was  noted  in  Reference  1 
that  the  predicted  additive  constant  is  about  2.5  compared  with  Patel's 
experimental  value  of  2.3. 

For  two-dimensional  flow  in  an  adverse  pressure  gradient  (Figures  3 
through  5)  there  is  a tendency  for  the  velocity  to  be  overpredicted 
in  the  inner  half  of  the  boundary  layer  --  the  shape  of  the  velocity 
profile  does  not  respond  sufficiently  to  the  effects  of  the  pressure 
gradient.  On  this  basis  one  would  presume  that  the  predicted  separa- 
tion point  would  lie  too  far  downstream.  On  the  other  hand,  for  at 
least  one  of  the  cases:  Schubauer  and  Spangenberg  (Figure  3),  the 
predicted  variation  of  with  x is  of  the  right  form  (even  though 
the  values  are  somewhat  low) . 

The  calculated  flow  over  the  swept  wing  was  matched  at  x/c  = 0.34 
(denoted  "measuring  station  #1"  in  Reference  14),  and  good  agreement 
between  theory  and  experiment  is  obtained  over  the  early  part  of  the 
run  (Figure  6) , except  for  the  underprediction  of  which  occurs  in 
most  of  the  comparisons.  At  x/c  = 0.54  (station  #4),  Figure  7 shows 
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that  the  predicted  vector  velocity  profile  is  in  excellent  agreement  with 
measurement,  and  that  the  predicted  profile  agrees  with  the  measured 
one  to  within  the  probable  level  of  experimental  accuracy.  Further  down- 
stream, however,  the  agreement  between  theory  and  experiment  is  less  good; 
in  fact  the  agreement  deteriorates  in  the  same  manner  as  other  investi- 
gators have  found  (see  Reference  16).  Figure  8,  which  corresponds  to 
x/c  = 0.74  (station  #7),  shows  the  familiar  overprediction  of  resultant 
velocity,  underprediction  of  crossflow  angle,  and  underprediction  of 
boundary- layer  thickness.  Interestingly,  the  turbulent-kinetic-energy 
profile  is  predicted  quite  well.  Separation,  defined  in  this  instance 
as  the  condition  where  the  limiting  streamlines  lie  parallel  to  the 
sweep  lines,  was  observed  to  occur  experimentally  at  about  x/c  = 0.87. 

The  calculation,  however,  continued  to  the  trailing  edge  without  any- 
where reaching  that  level  of  wall  streamline  deflection. 

The  calculation  for  the  supersonic  swept-shock  flow  was  adjusted  to 
give  the  appropriate  constant-pressure  profile  at  about  x = -l.in.,  where 
X is  measured  from  the  position  of  the  inviscid  shock  front.  At  this  con- 
dition, the  predicted  skin-friction  coefficient  lies  close  to  the  experi- 
mental value  given  by  a circular  Preston  tube.  The  latter  value  lies 
above  the  measurement  made  by  the  NAE  cobra  probe,  but  below  the  value 
derived  from  the  Spalding/Chi  correlation  [19]  for  two-dimensional  constant- 
pressure  boundary  layers  (Figure  9).  The  predicted  subsequent  decrease  of 
C^,  resulting  from  the  shock  interaction,  follows  a trend  which  is  consis- 
tent with  the  experimental  data.  However,  the  appreciable  scatter  of  the 
data  --  and  more  particularly  the  discrepancy  between  the  measurements 


9 


made  using  the  two  different  devices  --  precludes  any  detailed  comment 
about  the  accuracy  of  the  predictions.  The  wall  streamline  deflection 
angle  (measured,  here,  relative  to  the  tunnel  center  line)  is  shown  in 
the  lower  part  of  Figure  9.  The  prediction  of  this  quantity  was  found 
to  depend  rather  critically  on  the  assumed  conditions  along  the  edge  of 
the  integration  domain  nearest  to  the  shock  generator.  The  effect  of 
these  conditions,  on  the  solution  at  the  center  line,  conformed  with  the 

requirements  of  the  appropriate  zones  of  influence  but  was  still  appre-  j 

ciable.  The  solution  shown  in  Figure  9 corresponds  to  a conservative 

assumption  about  the  edge  conditions,  and  indicates  an  underprediction  i 

of  the  deflection  angle  on  the  center  line. 

The  predicted  crossflow  angle  through  the  boundary  layer  (measured 
relative  to  the  external  streamlines)  is  in  good  agreement  with  experi- 
ments except  close  to  the  wall  (Figure  10).  The  comparison  is  for  x = 0, 
which  is  three- fourths  the  way  up  the  pressure  rise.  The  resultant 
velocity  is  overpredicted  for  reasons  which  are  associated  with  the 
underprediction  of  streamline  deflection  angle.  The  predicted  total- 
temperature  profile  agrees  with  the  measured  one  to  within  2-3%;  the 
only  significant  discrepancy  appears  to  be  an  exaggeration  (in  the 
calculation)  of  the  variations  above  and  below  the  free-stream  total 
temperature.  The  predicted  and  measured  wall  temperatures  are  in  close 
agreement . 
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SECTION  IV.  MODIFICATIONS  TO  THE  MODEL 


Sensitivity  to  Changes  of  the  Empirical  Functions 

The  turbulence  model  involves  four  distinct  empirical  functions  which 
are  relevant  in  incompressible  flow; 

o the  diffusion  function  (Equation  2.5) 

o the  dissipation  length  (Equation  2.6) 

o the  Van  Driest  correction  to  dissipation  length  in 
the  wall  region  (Equation  2.9) 

o the  relationship  between  shear  stress  and  turbulent 
kinetic  energy  (Equations  2.11  and  2.12) 

It  has  previously  been  shown  [3]  that  the  form  of  the  diffusion  func- 
tion has  relatively  little  effect  on  the  solutions  and,  for  this  reason, 
it  was  left  unchanged  during  the  present  studies.  Effort  was  directed  to 
exploring  the  effects  of  modifying  the  other  empirical  functions.  Most 
of  these  sensitivity  studies  were  carried  out  using  incompressible  flat- 
plate  flow  as  a test  case.  Specifically,  the  following  modifications 
were  tried,  one  at  a time; 

(a)  changing  the  constant  in  the  numerator  of  the 
dissipation  length  function  (Equation  2.6)  to 
6.143,  which  reduces  the  values  of  through- 
out the  boundary  layer 

(b)  changing  the  denominator  of  the  dissipation 
length  function  to 

1 = Sn^  + I4n'^, 

which  reduces  the  values  of  over  the  central 
and  outer  parts  of  the  boundary  layer 

(c)  changing  the  denominator  of  the  dissipation 
length  function  to 

1 + 5n^  + n^ 

which  increases  the  values  of  Lp  over  the  outer 
part  of  the  boundary  layer 
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(d)  changing  the  value  of  A in  the  Van  Driest  formula  to 
13,  which  reduces  the  thickness  of  the  viscous  sub- 
layer 

(e)  increasing  the  constant  in  Equation  (2.13)  to  0.0278, 
which  increases  the  shear  stress  for  a given  level  of 
turbulent  intensity 

(f)  modifying  the  form  of  the  relationship  between  shear 
stress  and  turbulent  intensity  by  imposing  a ceiling 
of  0.15  on  the  value  of  f 

(g)  increasing  the  constant  in  Equation  (2.13)  to  0.0278 
and  also  imposing  a ceiling  of  0.167  on  the  value  of 
f. 

These  changes  are  illustrated  in  Figure  12.  In  modifications  (f)  and  (g) 
the  break  points  between  the  two  line  segments,  defining  the  function, 
were  arranged  to  lie  on  the  hyperbola  f = l/\p  which,  for  two-dimensional 
flow,  corresponds  to  the  statement  "production  = dissipation".  The  desir- 
ability of  imposing  the  ceiling  on  the  value  of  f was  suggested  partly  by 
experimental  data  and,  more  strongly,  by  conceptual  arguments  that  the 
value  of  f ought  to  reach  a maximum  for  large  rates  of  strain  --  or  even 
to  decrease  for  large  total  strains  [17].  A difficulty  arises,  however, 
that  any  departure  from  a proportionality  relationship  between  f and  ip 
implies  that  the  resulting  model  is  not  invariant  to  rotation  of  the 
coordinate  axes.  There  are  ways  of  surmounting  this  difficulty,  by 
redefining  the  parameters  involved,  but  they  were  not  considered  here 
because  significant  coding  changes  would  have  been  required.  The  results 
obtained  from  modifications  (f)  and  (g) , above,  therefore  have  to  be 
evaluated  in  terms  of  the  particular  coordinate  system  used  here  (parallel 
and  normal  to  the  undistrubed  flow  at  infinity) . 
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The  results  of  the  sensitivity  study  are  shovm  in  Figures  12  through  14. 
Here  the  predicted  velocity  profiles,  and  profiles  of  q^,  are  compared  for  a 
Reynolds  number  of  10**  based  on  displacement  thickness.  The  main  observation 
is  that  none  of  the  changes  has  more  than  a minimal  effect  on  the  velocity 
profile,  but  some  of  them  have  a pronounced  effect  on  turbulent  kinetic 
energy.  The  effect  on  turbulent  kinetic  energy  is  surprisingly  large  in 
the  case  of  the  modification  to  the  Van  Driest  constant.  It  might  be 
thought  that  this  constant  would  only  affect  conditions  in  the  sublayer; 
however,  it  emerges  that  the  effect  is  distributed  across  the  whole  boundary 
layer.  The  interpretation  is  that  the  level  of  at  the  edge  of  the  sub- 
layer through  which  the  variation  of  q^  is  rapid  serves  as  a boundary 
condition  for  the  solution  throughout  the  fully  turbulent  region  where 
3q^/3y  is  relatively  small. 

In  addition  to  the  results  presented  in  Figures  12  through  14,  the 
following  effects  of  the  modifications  on  wall  shear  stress  were  observed 

(these  are  all  for  a value  of  R = 10**) : 

®6* 

Modification  *^f 


(a) 

0.00212 

(b) 

0.00244 

(c) 

0.00248 

(d) 

0.00296 

(e) 

0.00298 

(f) 

0.00234 

(g) 

0.00270 

These  values  may  be  compared  with  0.00256,  obtained  using  the  original 
model,  and  with  0.0028  which  is  Klebanoff's  measured  value. 
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The  changes  to  the  empirical  functions,  considered  so  far,  have  been 
evaluated  in  terms  of  their  effects  in  incompressible  flow.  No  separate 
sensitivity  study  has  been  conducted  under  conditions  of  compressibility. 

One  of  the  underlying  assumptions  is  that  the  turbulence  model  remains 
valid  until  the  Mach  number  based  on  (q^)^  becomes  comparable  with  unity, 
and  that  is  unlikely  to  occur  below  a free-stream  Mach  number  of  about  3. 

The  only  empirical  function  that  must  be  studied  in  a context  of  compress- 
ible flow  is  the  relationship:  between  turbulent  heat  flux  and 

turbulent  kinetic  energy  (Equation  2.14).  A few  calculations  were  done  in 
which  a ceiling  was  imposed  on  the  value  of  g,  by  analogy  with  modifications 
(f)  and  (g) , above.  The  results  indicated  that  such  a change  was  inadmis- 
sible because  it  destroyed  the  conservation  of  total  enthalpy  in  the 
boundary  layer:  the  average  total  temperature  was  found  to  exceed  the 
free-stream  value.  The  only  modification,  made  to  the  relationship  in  the 
present  work,  was  to  change  the  constant  in  the  numerator  of  Equation  (2.15): 
i.e.  0.0225,  to  0.0278  whenever  the  latter  value  was  used  in  Equation  (2.13). 
In  this  way  the  assumed  value  of  the  turbulent  Prandtl  number  was  preserved, 
at  least  for  small  values  of  i|<,  regardless  of  the  form  chosen  for  f(ij;). 
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History  Effect  on  Dissipation  Length 


In  addition  to  the  modifications  considered  above,  some  studies  were 
made  to  explore  the  effects  of  uncoupling  the  dissipation  length  from  the 
local  boundary- layer  thickness.  There  are  strong  conceptual  arguments,  and 
some  experimental  evidence  suggesting  that  the  scale  of  the  large  eddies, 
which  determines  the  dissipation  length,  responds  relatively  slowly  to 
changes  of  boundary- layer  thickness.  The  assumption,  in  the  original  model, 
that  this  response  is  immediate,  leads  to  an  overestimate  of  dissipation 
length  in  boundary  layers  (such  as  that  of  van  den  Berg,  et.  al . [14,  15]) 
which  are  characterized  by  large  spatial  rates  of  growth. 


In  order  to  reflect  the  finite  rate  of  response  of  the  dissipation 
length,  it  is  necessary  to  incorporate  a second  turbulence  equation.  Some 
suggestions  are  given  in  Reference  18  as  to  the  appropriate  form  of  such 
an  equation*,  but  here  the  global  effect  is  represented  simply  by  replacing 
6,  in  Equation  (2.6),  by  6^  and  by  redefining  n (in  the  same  equation)  as 


n 


(16) 


whers  6^  is  an  effective  thickness  which  lags  behind  the  physical  thickness, 
6,  according  to  the  rate  equation 


96, 


96, 


U + W (6  - 6 ) . 

e 9x  e 9z  k6  e 


(15) 


The  rate  constant,  k,  is  taken  provisionally  to  be  6,  implying  that  the 


Other  investigators  have  proposed  appropriate  forms  of  an 
equation  for  dissipation  rate. 
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dissipation  length  in  the  outer  part  of  the  layer  has  a characteristic 
response  time  equal  to  the  time  taken  for  a fluid  particle  to  travel 
a distance  of  66.  In  the  inner  part  of  the  boundary  layer  the  charac- 
teristic time  becomes  shorter,  and  it  reaches  zero  at  the  wall.  The 
form  of  Equation  (2.6),  with  n now  given  by  Equation  (4.1),  ensures 
a qualitatively  correct  variation  of  the  characteristic  time  across 
the  boundary  layer. 


» 
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Proposed  Changes  to  the  Turbulence  Model 

The  most  significant  weakness  of  the  original  model  appears  to  lie  in 
its  lack  of  responsiveness  to  adverse  pressure  gradients  leading  to  separa- 
tion. The  simplest  means  of  correcting  this  weakness  is  to  modify  the  rela- 
tionship between  shear  stress  and  turbulent  kinetic  energy.  Imposing  a 
ceiling  on  the  value  of  f (Equations  2.11  and  2.12)  prevents  the  generation 
of  the  high  shear-stress  levels  which,  in  large  measure,  cause  the  unres- 
ponsiveness. Modification  (f) , defined  above  provides  a satisfactory  1 

provisional  improvement  although,  as  noted  earlier,  it  suffers  from  lack 
of  invariance  to  rotation  of  the  coordinate  system.  It  also  suffers  from 
the  deficiency  of  producing  wall  shear-stress  values  which  typically  are 
too  low.  For  the  flat-plate  flow  (Section  III,  above)  it  was  noted  that 
C_  lay  about  16%  below  the  experimental  value,  for  R = 10**.  This 
modification  will  be  referred  to  as  "Proposed  Model  A"  in  what  follows. 

An  attempt  has  been  made  to  derive  a second  model  which  would  exhibit 
the  proper  degree  of  responsiveness,  while  also  producing  the  correct  level 
of  wall  shear  stress.  Considerably  more  effort  was  put  into  the  derivation 
of  this  model:  denoted  "Proposed  Model  B",  and  it  seems  to  give  satisfactory 
agreement  with  all  the  cases  considered.  Model  B reflects  the  following 
combination  of  changes  to  the  empirical  functions: 

o Modification  of  the  dissipation  length 

function  (modification  (c)  of  Section 
IV,  above) 

o Modification  of  the  shear-stress/turbulent- 

kinetic-energy  relationship  (modification 
(g)  of  Section  IV) 

o Incorporation  of  the  history  effect  on  dissi- 
pation length  via  Equation  (4.2) 
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The  proposed  modification  to  the  dissipation  length  (Figure  11)  has 
the  effect  of  increasing  the  levels  of  and  shear  stress  in  the  outer 
part  of  the  boundary  layer.  This  is  desirable,  for  flows  in  small  adverse 
pressure  gradients,  to  offset  the  decrease  brought  about  by  allowing  the 
dissipation  length  to  lag  behind  the  local  boundary- layer  thickness.  For 
flows  in  severe  adverse  pressure  gradients  the  two  effect  do  not  offset 
one  another  and  the  proper  reduction  in  shear  stress  is  achieved.  The 
modification  to  the  f-function  (Figure  11)  further  limits  the  growth  of 
turbulent  kinetic  energy  in  strong  adverse  pressure  gradients.  However, 
the  higher  slope  of  the  function,  for  small  values  of  the  argument,  im- 
proves the  flat-plate  solution  and  results  in  more  accurate  wall  shear- 
stress  predictions.  It  was  not  found  desirable  to  change  the  Van  Driest 
constant;  nor,  as  mentioned  earlier,  have  any  possible  modifications  to 
the  diffusion  function  been  explored. 
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Comparison  Between  the  Proposed  Models  and  Experiment 

Calculations  involving  the  Proposed  Models  A and  B are  compared  with 
the  original  calculations,  and  with  the  experimental  data,  in  Figures  15 
through  23. 


For  two-dimensional  flow  over  a flat  plate  there  was  little  room  for 
improvement  as  far  as  the  profiles  of  velocity  and  are  concerned. 
Figure  15  shows  that  agreement  with  the  measured  velocity  profile  may 


perhaps  be  slightly  better  with  the  new  models,  whereas  agreement  with 
is  slightly  worse  (if  the  measurements  are  sufficiently  accurate  to 


make  such  a judgment).  The  associated  wall  shear-stress  values,  for  a 


Reynolds  number  of  10**  based  on  displacement  thickness,  are  as  follows: 


Model 

Original  Model 
Proposed  Model  A 
Proposed  Model  B 
Experiment  (Klebanoff) 
Experiment  (Smith  § Walker) 


0.00256 

0.00234 

0.00272 

0.0028 

0.0027 


There  is  some  evidence  that  the  error  in  C^  for  Model  A results 


partly  from  truncation  errors  in  the  solution  method,  and  is  substantially 


reduced  when  the  calculation  is  done  to  second-order  accuracy.  It  should 
be  recognized  that  the  fine-tuning  of  the  turbulence  model,  attempted  here, 
may  have  been  hampered  by  numerical  errors  (although  these  are  believed 
generally  to  be  small).  Furthermore,  the  possibility  exists  that  adjust- 
ments to  the  turbulence  model  might  have  been  made  simply  to  offset  errors 
in  the  numerical  scheme,  and  that  these  adjustments  will  prove  to  have  been 
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unnecessary  when  verification  of  the  second-order  version  of  the  method  is 
completed. 


The  two  new  models  appear  to  have  corrected  the  lack  of  responsiveness 
to  adverse  pressure  gradients  leading  to  separation.  Agreement  with  the  two 
sets  of  two-dimensional  data  (Figures  16  through  18)  is  improved,  and  the 
agreement  with  the  swept-wing  data  is  substantially  improved  (Figures  19 
through  21) . Model  B,  in  particular,  produces  approximately  the  correct 
variation  of  crossflow  angle  through  the  boundary  layer,  at  Station  #7. 

Both  new  models  produce  good  agreement  with  the  measured  resultant  velocity 
profile,  although  there  is  a tendency  to  underpredict  the  level  of  turbulent 
kinetic  energy.  As  in  the  case  of  the  flat -plate  flow.  Model  B appears  to 
give  the  most  accurate  predictions  of  wall  shear  stress. 

The  proposed  modifications  have  only  a minimal  effect  on  the  comparison 
with  the  swept-shock  data  (Figures  22  and  23) . This  observation  seems  to 
indicate  that  the  discrepancies  between  the  original  model  and  Peake’s 
experiment  were  not  closely  related  to  the  discrepancies  between  it  and 
the  other  test  cases.  The  reasons  behind  the  unsatisfactory  agreement  with 
the  swept-shock  data  need  to  be  explored  more  fully,  and  a more  sophisti- 
cated calculation  performed  in  which  edge  boundary  conditions  can  be  con- 
trolled more  precisely.  As  it  stands  the  comparison  does  not  provide  either 
positive  or  negative  support  for  the  proposed  changes  to  the  model. 


I 
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SECTION  V.  CONCLUSIONS  AND  RECOMMENDATIONS 


Subject  to  rather  severe  constraints  imposed  on  this  study  by  funding 
limi.  .ons,  the  immediate  objectives  have  been  accomplished.  The  turbu- 
lence model,  on  which  the  method  of  Reference  1 is  based,  has  been  found  to 
behave  well  except  close  to  separation.  In  small  and  moderate  pressure 
gradients  the  predicted  profiles  of  velocity  and  turbulent  kinetic  energy 
are  in  satisfactory  agreement  with  experiment.  At  low  Reynolds  numbers  and 
Mach  numbers  the  predicted  wall  shear  stress  is  within  about  5 - 10%  of 
measured  values,  and  the  error  appears  to  decrease  with  increase  of  either 
Reynolds  number  or  Mach  number.  On  the  other  hand,  the  model  exhibits  a 
lack  of  responsiveness  to  strong  adverse  pressure  gradients  leading  to 
separation,  and  this  weakness  seems  to  be  somewhat  more  serious  for  three- 
dimensional  flows  than  for  those  in  two  dimensions.  For  one  of  the  three- 
dimensional  flows  studied  here:  flow  over  a swept  wing,  the  streamline 
deflections  were  underpredicted,  and  correspondingly,  the  observed  pattern 
of  separation  onset  was  not  reflected  by  the  calculation.  For  the  rather 
extreme  case  of  supersonic  flow  through  a swept  shock,  the  model  performed 
as  well  as  could,  perhaps,  be  expected.  Discrepancies  with  the  data  were 
found  in  this  case;  however  the  importance  of  the  discrepancies  is  difficult 
to  assess  because  of  uncertainties  in  the  measurements  and  the  difficulty 
of  setting  up  the  calculation  to  reflect  appropriate  boundary  conditions. 

Attempts  have  been  made  to  fine-tune  the  turbulence  model  by  making  a 
number  of  modifications  and  observing  the  consequent  effect  on  the  agree- 
ment between  theory  and  experiment.  The  results  are  not  conclusive,  primarily 
because  the  optimization  involved  only  a few  sets  of  data,  but  two  alternative 
new  models  (represented  by  particular  choices  of  the  empirical  functions)  are 
proposed  which  seem  to  perform  better  than  the  original  model. 
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One  of  them  consists  of  a modification  to  the  functional  relationship 
between  shear  stress  and  turbulent  kinetic  energy.  The  other  consists  of 
a modification  to  this  same  relationship  and  also  a modification  to  the  way 
in  which  the  dissipation  length  is  represented.  In  the  second  model  the 
dissipation  length  is  expressed,  as  before,  as  a function  of  position  through 
the  boundary  layer,  but  now  relative  to  a characteristic  length  which  lags 
behind  the  local  boundary- layer  thickness.  This  lag  is  intended  to  reflect 
the  finite  response  time  governing  the  dynamic  relationship  between  the  size 
of  the  large  eddies  and  the  scale  of  the  turbulent  flow. 

Implementation  of  the  second  of  these  proposed  modifications  in  the 
computer  code  has  only  been  carried  out  rather  crudely,  and  both  modifi- 
cations involve  a weakness  regarding  invariance  to  the  coordinate  system. 

The  shear-stress/turbulent-kinetic-energy  relationship  corresponds  to  a 
generalized  constitutive  relationship  between  stress  and  strain.  The  form 
provisionally  adopted  in  the  proposed  modifications  corresponds  to  a con- 
stitutive relationship  which  does  not  have  the  necessary  property  of  in- 
variance. A further  modification  would  be  needed  to  recover  this  property, 
and  certain  associated  coding  changes  would  be  involved.  It  is  strongly 
recommended  that  this  additional  work  be  carried  out  so  as  to  improve  the 
generality  of  the  model. 

It  is  worth  mentioning,  too,  that  the  study  of  alternate  constitutive 
relationships  would  likely  be  a fruitful  area  of  research  in  turbulent 
flows.  The  limitations  of  isotropic  eddy- viscosity  models  (even  those, 
like  the  model  of  Reference  1,  in  which  the  magnitude  of  the  eddy  viscosity 
is  not  prescribed  in  advance)  are  well  known.  But  little  has  been  done  to 
remedy  their  deficiencies  other  than  certain  attempts  to  calculate  all  the 
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elements  of  the  Reynolds  stress  tensor  separately.  Such  an  approach 
introduces  even  more  uncertainties,  and  it  may  be  preferable,  insteady,  to 
represent  the  stress  tensor  in  terms  of  a smaller  number  of  independent 
variables  together  with  an  appropriate  constitutive  relationship  for 
generating  the  individual  elements.  It  is  rapidly  becoming  apparent  that 
the  normal  stresses  can  not  be  ignored  if  the  model  is  to  have  general 
validity,  and  a sound,  economical  means  of  dealing  with  the  complete 
tensor  is  urgently  needed. 

Another  urgent  extension  of  this  study  would  consist  of  testing  the 
proposed  modifications  to  the  model  in  a second-order-accurate  version  of 
the  calculation  method.  This  step  is  necessary  to  check  whether  the  modi- 
fications actually  represent  improvements  in  the  turbulence  model,  or 
represent  nothing  more  than  devices  for  generating  additional  errors  to 
offset  those  associated  with  the  numerical  scheme. 
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X = 3.132  (meters) 
(match  at  X = 1.282  m.) 


FIGURE  4 COMPARISON  WITH  THE  DATA  OF  LUDWIEG  AND 

TILLMANN  (VELOCITY  PROFILE) 
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FIGURE  5 


COMPARISON  WITH  THE  DATA  OF  SCHUBAUER  AND 
SPANGENBERG  (VELOCITY  PROFILE) 


FIGURE  10  COMPARISON  WITH  THE  DATA  OF  PEAKE 

(VELOCITY  AND  TEMPERATURE  PROFILES) 
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original  model 


FIGURE  16 


COMPARISON  WITH  THE  DATA  OF  LUDWIEG  AND  TILLMANN 
AND  SCHUBAUER  AND  SPANGENBERG  (WALL  SHEAR  STRESS 
AND  BOUNDARY- LAYER  THICKNESS) 
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FIGURE  17 


COMPARISON  WITH  THE  DATA  OF  LUDWIEG  AND  TILLMANN 
(VELOCITY  PROFILE) 


FIGURE  18 


COMPARISON  WITH  THE  DATA  OF  SCHUBAUER  AND  SPANGENBERG 
(VELOCITY  PROFILE) 
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FIGURE  19  COMPARISON  WITH  THE  DATA  OF  VAN  DEN  BERG  AND 

ELSENAAR  (WALL  SHEAR  STRESS  AND  BOUNDARY -LAYER 
THICKNESS) 
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FIGURE  23  COMPARISON  WITH  THE  DATA  OF  PEAKE 

(VELOCITY  AND  TEMPERATURE  PROFILES) 
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